1/f Noise Project Manuscript 1 Final Analysis
Fine-Scale Spatial Heterogeneity in Phytochemistry Increases Herbivore Performance by Enhancing Herbivore Foraging
Vincent S. Pan, Enakshi Ghosh, Paul J. Ode, William C. Wetzel, Kadeem J. Gilbert, Ian S. Pearse
May 23, 2024
1 Setup
1.1 Load libraries
library(survival)
library(piecewiseSEM)
source("spat1f/init_analysis.R")Data last updated on 2024-05-031.2 Handy functions not in spat1f
trend_by_pre_wt <- function(model, term){
emmeans::emtrends(model, var = "cat_pre_wt_log_scale", specs = term) %>%
pairs()
}
contrast_by_pre_wt <- function(model, term, type = c("pairwise", "trt.vs.ctrl"), at = c(-2, 2)){
type <- match.arg(type)
f <- as.formula(sprintf("%s ~ %s | cat_pre_wt_log_scale", type, term))
emmeans::emmeans(model,
f,
var = term,
at = list(cat_pre_wt_log_scale = at))
}
check_mod <- function(model){
DHARMa::simulateResiduals(model) %>%
plot()
performance::posterior_predictive_check(model)
}1.3 Quick cleaning
# Main data.frame for analyses
d <- ref_data %>%
filter(error == 0) %>%
mutate(
cat_pre_wt_log = log(cat_pre_wt),
cat_pre_wt_log_scale = as.numeric(scale(log(cat_pre_wt))),
cat_size = ifelse(cat_pre_wt < median(cat_pre_wt), "small", "big"),
has_var = ifelse(var_trt != "constant", 1, 0),
mean_toxic_conc_scale = as.numeric(scale(mean_toxic_conc)),
mean_toxic_scale = as.numeric(scale(mean_toxic)),
area_herb_log_scale = as.numeric(scale(log(area_herb+1))),
var_toxic_12_scale = as.numeric(scale(var_toxic_12)),
on_toxic_logit_scale = as.numeric(scale(adjust_prop(on_toxic,
trans = "emp",
nudge.method = "none",
na.action = "ignore"))),
ava_qual_scale = as.numeric(scale(ava_qual)),
ava_qual_logit_scale = as.numeric(scale(adjust_prop(ava_qual,
trans = "emp",
nudge.method = "none",
na.action = "ignore"))),
sl_mean_obs_log_scale = as.numeric(scale(log(sl_mean_obs))),
cat_pre_wt_log_scale_sq = as.numeric(scale(log(cat_pre_wt)^2)),
RGR_scale = as.numeric(scale(RGR)),
sl_mean_pred1 = shape_1 * scale_1,
sl_mean_pred2 = shape_2 * scale_2,
sl_mean_pred3 = shape_3 * scale_3,
sl_mean_pred4 = shape_4 * scale_4,
k1 = scale_1 / scale_2,
k2 = scale_3 / scale_4,
prop_explore_logit = qlogis(prop_explore),
prop_explore_logit_scale = as.numeric(scale(prop_explore_logit))
)
# Subset of data.frame for SEM
d2 <- d %>%
filter(
var_trt != "constant"
) %>%
filter(
!is.na(area_herb_log_scale) &
!is.na(var_toxic_12_scale) &
!is.na(mean_toxic_conc_scale) &
!is.na(sl_mean_obs) &
!is.na(prop_explore) &
!is.na(RGR)) %>%
mutate(
var_high = ifelse(var_trt == "high_var", 1, 0),
beta_red = ifelse(as.character(beta) == 5, 1, 0),
beta_white = ifelse(as.character(beta) == 0, 1, 0),
beta_red_cat_pre_wt_log_scale = beta_red * cat_pre_wt_log_scale,
beta_white_cat_pre_wt_log_scale = beta_white * cat_pre_wt_log_scale,
var_high_cat_pre_wt_log_scale = var_high * cat_pre_wt_log_scale,
beta_numeric_scale = as.numeric(scale(as.numeric(beta)))
)
# Subset of data.frame for more detailed analyses
d3 <- d %>%
filter(
var_trt != "constant"
) %>%
filter(
sl_mean_pred2 < 2000 &
sl_mean_pred1 > 0 &
sl_mean_pred2 > 0 &
sl_mean_pred3 > 0 &
sl_mean_pred4 > 0 &
scale_1 > 0 &
shape_1 > 0 &
scale_2 > 0 &
shape_2 > 0 &
scale_3 > 0 &
shape_3 > 0 &
scale_4 > 0 &
shape_4 > 0 &
!is.na(area_herb_log_scale) &
!is.na(var_toxic_12_scale) &
!is.na(mean_toxic_conc_scale) &
!is.na(on_toxic_logit_scale) &
!is.na(ava_qual_scale) &
!is.na(sl_mean_obs) &
!is.na(prop_explore) &
!is.na(RGR)) %>%
mutate(
var_high = ifelse(var_trt == "high_var", 1, 0),
beta_red = ifelse(as.character(beta) == 5, 1, 0),
beta_white = ifelse(as.character(beta) == 0, 1, 0),
beta_red_cat_pre_wt_log_scale = beta_red * cat_pre_wt_log_scale,
beta_white_cat_pre_wt_log_scale = beta_white * cat_pre_wt_log_scale,
var_high_cat_pre_wt_log_scale = var_high * cat_pre_wt_log_scale,
beta_numeric_scale = as.numeric(scale(as.numeric(beta)))
)2 Herbivore Performance
2.1 RGR
2.1.1 Model summary
rgr_m <- glmmTMB(
RGR ~
(var_trt + beta) * cat_pre_wt_log_scale +
I(cat_pre_wt_log_scale^2) + # residual plot show significant non-linearity
(1|session_id),
data = d %>%
filter(
var_trt != "constant"
)
); summary(rgr_m) Family: gaussian ( identity )
Formula:
RGR ~ (var_trt + beta) * cat_pre_wt_log_scale + I(cat_pre_wt_log_scale^2) +
(1 | session_id)
Data: d %>% filter(var_trt != "constant")
AIC BIC logLik deviance df.resid
-879.2 -849.4 450.6 -901.2 100
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
session_id (Intercept) 3.979e-06 0.001995
Residual 1.603e-05 0.004004
Number of obs: 111, groups: session_id, 5
Dispersion estimate for gaussian family (sigma^2): 1.6e-05
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 0.0162504 0.0013457 12.075 < 2e-16 ***
var_trtlow_var -0.0034872 0.0007740 -4.506 6.62e-06 ***
beta0 -0.0003065 0.0009664 -0.317 0.75112
beta5 -0.0013545 0.0009173 -1.477 0.13978
cat_pre_wt_log_scale -0.0048759 0.0011991 -4.066 4.78e-05 ***
I(cat_pre_wt_log_scale^2) -0.0021953 0.0006731 -3.262 0.00111 **
var_trtlow_var:cat_pre_wt_log_scale 0.0025873 0.0007897 3.276 0.00105 **
beta0:cat_pre_wt_log_scale 0.0011924 0.0009491 1.256 0.20897
beta5:cat_pre_wt_log_scale 0.0024731 0.0009710 2.547 0.01087 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1car::Anova(rgr_m, type = "III")Analysis of Deviance Table (Type III Wald chisquare tests)
Response: RGR
Chisq Df Pr(>Chisq)
(Intercept) 145.8141 1 < 2.2e-16 ***
var_trt 20.3008 1 6.617e-06 ***
beta 2.3963 2 0.301758
cat_pre_wt_log_scale 16.5348 1 4.777e-05 ***
I(cat_pre_wt_log_scale^2) 10.6374 1 0.001108 **
var_trt:cat_pre_wt_log_scale 10.7334 1 0.001052 **
beta:cat_pre_wt_log_scale 6.4936 2 0.038899 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 12.1.2 Post-hoc contrasts
trend_by_pre_wt(rgr_m, "var_trt") contrast estimate SE df t.ratio p.value
high_var - low_var -0.00259 0.00079 100 -3.276 0.0014
Results are averaged over the levels of: beta contrast_by_pre_wt(rgr_m, "var_trt")$emmeans
cat_pre_wt_log_scale = -2:
var_trt emmean SE df lower.CL upper.CL
high_var 0.014224 0.00355 100 0.007178 0.02127
low_var 0.005562 0.00312 100 -0.000635 0.01176
cat_pre_wt_log_scale = 2:
var_trt emmean SE df lower.CL upper.CL
high_var -0.000393 0.00208 100 -0.004512 0.00373
low_var 0.001295 0.00238 100 -0.003421 0.00601
Results are averaged over the levels of: beta
Confidence level used: 0.95
$contrasts
cat_pre_wt_log_scale = -2:
contrast estimate SE df t.ratio p.value
high_var - low_var 0.00866 0.00181 100 4.790 <.0001
cat_pre_wt_log_scale = 2:
contrast estimate SE df t.ratio p.value
high_var - low_var -0.00169 0.00171 100 -0.988 0.3256
Results are averaged over the levels of: beta trend_by_pre_wt(rgr_m, "beta") contrast estimate SE df t.ratio p.value
(beta-5) - beta0 -0.00119 0.000949 100 -1.256 0.4232
(beta-5) - beta5 -0.00247 0.000971 100 -2.547 0.0329
beta0 - beta5 -0.00128 0.000940 100 -1.363 0.3644
Results are averaged over the levels of: var_trt
P value adjustment: tukey method for comparing a family of 3 estimates contrast_by_pre_wt(rgr_m, "beta")$emmeans
cat_pre_wt_log_scale = -2:
beta emmean SE df lower.CL upper.CL
-5 0.012890 0.00382 100 0.005313 0.02047
0 0.010199 0.00333 100 0.003591 0.01681
5 0.006589 0.00319 100 0.000255 0.01292
cat_pre_wt_log_scale = 2:
beta emmean SE df lower.CL upper.CL
-5 -0.001439 0.00226 100 -0.005922 0.00304
0 0.000639 0.00231 100 -0.003940 0.00522
5 0.002153 0.00257 100 -0.002946 0.00725
Results are averaged over the levels of: var_trt
Confidence level used: 0.95
$contrasts
cat_pre_wt_log_scale = -2:
contrast estimate SE df t.ratio p.value
(beta-5) - beta0 0.00269 0.00221 100 1.218 0.4456
(beta-5) - beta5 0.00630 0.00223 100 2.820 0.0158
beta0 - beta5 0.00361 0.00211 100 1.711 0.2062
cat_pre_wt_log_scale = 2:
contrast estimate SE df t.ratio p.value
(beta-5) - beta0 -0.00208 0.00205 100 -1.015 0.5689
(beta-5) - beta5 -0.00359 0.00206 100 -1.746 0.1935
beta0 - beta5 -0.00151 0.00211 100 -0.719 0.7529
Results are averaged over the levels of: var_trt
P value adjustment: tukey method for comparing a family of 3 estimates 2.1.3 Misc stats
r2(rgr_m, tolerance = 10^-10)# R2 for Mixed Models
Conditional R2: 0.555
Marginal R2: 0.4442.1.4 Check model
check_mod(rgr_m)
2.1.5 Plot
g1 <- marginal_effects(rgr_m, terms = c("cat_pre_wt_log_scale","beta")) %>%
ggplot(aes(x = unscalelog(d$cat_pre_wt_log)(cat_pre_wt_log_scale), y = yhat)) +
geom_ribbon(aes(ymax = upper, ymin = lower, fill = beta), alpha = 0.2) +
geom_line(aes(color = beta), linewidth = 2) +
geom_point(
data = rgr_m$frame,
aes(x = unscalelog(d$cat_pre_wt_log)(cat_pre_wt_log_scale), y = RGR, color = beta),
size = 3, alpha = 0.8
) +
theme_bw(base_size = 15) +
theme(legend.position = "top") +
scale_x_continuous(trans = "log10", labels = fancy_scientific) +
scale_color_brewer(type = "qual", aesthetics = c("fill","color")) +
labs(x = "Cat pre-weight (g)", y = expression(RGR~(hour^-1)),
color = expression(beta), fill = expression(beta))
g2 <- marginal_effects(rgr_m, terms = c("cat_pre_wt_log_scale","var_trt")) %>%
ggplot(aes(x = unscalelog(d$cat_pre_wt_log)(cat_pre_wt_log_scale), y = yhat)) +
geom_ribbon(aes(ymax = upper, ymin = lower, fill = var_trt), alpha = 0.2) +
geom_line(aes(color = var_trt), linewidth = 2) +
geom_point(
data = rgr_m$frame,
aes(x = unscalelog(d$cat_pre_wt_log)(cat_pre_wt_log_scale), y = RGR, color = var_trt),
size = 3, alpha = 0.8
) +
theme_bw(base_size = 15) +
scale_x_continuous(trans = "log10", labels = fancy_scientific) +
theme(legend.position = "top") +
scale_color_discrete(type = c("navy","skyblue"),
aesthetics = c("fill","color"),
label = c("High", "Low")) +
labs(x = "Cat pre-weight (g)", y = expression(RGR~(hour^-1)),
color = "Variation", fill = "Variation")
ggarrange(g1, g2 + labs(y = ""))
2.2 Time to pupation
2.2.1 Model summary
pupt_m_full <- glmmTMB(
time_to_pupation ~
(var_trt + beta) * cat_pre_wt_log_scale +
I(cat_pre_wt_log_scale^2) + # Residuals show significant non-linearity
(1|session_id),
data = d %>%
filter(var_trt != "constant"),
family = Gamma(link = "log"),
); summary(pupt_m_full) Family: Gamma ( log )
Formula: time_to_pupation ~ (var_trt + beta) * cat_pre_wt_log_scale +
I(cat_pre_wt_log_scale^2) + (1 | session_id)
Data: d %>% filter(var_trt != "constant")
AIC BIC logLik deviance df.resid
338.1 366.3 -158.0 316.1 85
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
session_id (Intercept) 2.839e-11 5.328e-06
Number of obs: 96, groups: session_id, 5
Dispersion estimate for Gamma family (sigma^2): 0.0254
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 2.180302 0.036853 59.16 < 2e-16 ***
var_trtlow_var 0.102145 0.033998 3.00 0.00266 **
beta0 -0.001845 0.040823 -0.05 0.96395
beta5 0.059543 0.041181 1.45 0.14821
cat_pre_wt_log_scale -0.428662 0.032643 -13.13 < 2e-16 ***
I(cat_pre_wt_log_scale^2) -0.081267 0.019069 -4.26 2.03e-05 ***
var_trtlow_var:cat_pre_wt_log_scale 0.026518 0.038876 0.68 0.49516
beta0:cat_pre_wt_log_scale -0.030363 0.043195 -0.70 0.48210
beta5:cat_pre_wt_log_scale -0.126801 0.047095 -2.69 0.00709 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1car::Anova(pupt_m_full, type = "III")Analysis of Deviance Table (Type III Wald chisquare tests)
Response: time_to_pupation
Chisq Df Pr(>Chisq)
(Intercept) 3500.2207 1 < 2.2e-16 ***
var_trt 9.0264 1 0.002661 **
beta 2.7351 2 0.254736
cat_pre_wt_log_scale 172.4448 1 < 2.2e-16 ***
I(cat_pre_wt_log_scale^2) 18.1614 1 2.029e-05 ***
var_trt:cat_pre_wt_log_scale 0.4653 1 0.495162
beta:cat_pre_wt_log_scale 7.5589 2 0.022836 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1pupt_m <- glmmTMB(
time_to_pupation ~
var_trt + beta * cat_pre_wt_log_scale +
I(cat_pre_wt_log_scale^2) + # Residuals show significant non-linearity
(1|session_id),
data = d %>%
filter(var_trt != "constant"),
family = Gamma(link = "log"),
); summary(pupt_m) Family: Gamma ( log )
Formula:
time_to_pupation ~ var_trt + beta * cat_pre_wt_log_scale + I(cat_pre_wt_log_scale^2) +
(1 | session_id)
Data: d %>% filter(var_trt != "constant")
AIC BIC logLik deviance df.resid
336.5 362.2 -158.3 316.5 86
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
session_id (Intercept) 2.159e-11 4.647e-06
Number of obs: 96, groups: session_id, 5
Dispersion estimate for Gamma family (sigma^2): 0.0255
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 2.180555 0.036937 59.03 < 2e-16 ***
var_trtlow_var 0.107447 0.033197 3.24 0.00121 **
beta0 -0.003578 0.040842 -0.09 0.93019
beta5 0.059536 0.041271 1.44 0.14915
cat_pre_wt_log_scale -0.420147 0.030221 -13.90 < 2e-16 ***
I(cat_pre_wt_log_scale^2) -0.083365 0.018881 -4.42 1.01e-05 ***
beta0:cat_pre_wt_log_scale -0.026011 0.042849 -0.61 0.54383
beta5:cat_pre_wt_log_scale -0.126834 0.047183 -2.69 0.00719 **
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1car::Anova(pupt_m, type = "III")Analysis of Deviance Table (Type III Wald chisquare tests)
Response: time_to_pupation
Chisq Df Pr(>Chisq)
(Intercept) 3485.0791 1 < 2.2e-16 ***
var_trt 10.4759 1 0.001209 **
beta 2.8053 2 0.245947
cat_pre_wt_log_scale 193.2743 1 < 2.2e-16 ***
I(cat_pre_wt_log_scale^2) 19.4946 1 1.009e-05 ***
beta:cat_pre_wt_log_scale 7.6890 2 0.021397 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 12.2.2 Post-hoc contrasts
contrast_by_pre_wt(pupt_m, "var_trt")$emmeans
cat_pre_wt_log_scale = -2:
var_trt emmean SE df asymp.LCL asymp.UCL
high_var 2.808 0.0769 Inf 2.657 2.96
low_var 2.915 0.0818 Inf 2.755 3.08
cat_pre_wt_log_scale = 2:
var_trt emmean SE df asymp.LCL asymp.UCL
high_var 0.924 0.0672 Inf 0.792 1.06
low_var 1.031 0.0725 Inf 0.889 1.17
Results are averaged over the levels of: beta
Results are given on the log (not the response) scale.
Confidence level used: 0.95
$contrasts
cat_pre_wt_log_scale = -2:
contrast estimate SE df z.ratio p.value
high_var - low_var -0.107 0.0332 Inf -3.237 0.0012
cat_pre_wt_log_scale = 2:
contrast estimate SE df z.ratio p.value
high_var - low_var -0.107 0.0332 Inf -3.237 0.0012
Results are averaged over the levels of: beta
Results are given on the log (not the response) scale. trend_by_pre_wt(pupt_m, "beta") contrast estimate SE df z.ratio p.value
(beta-5) - beta0 0.026 0.0428 Inf 0.607 0.8163
(beta-5) - beta5 0.127 0.0472 Inf 2.688 0.0197
beta0 - beta5 0.101 0.0473 Inf 2.132 0.0834
Results are averaged over the levels of: var_trt
P value adjustment: tukey method for comparing a family of 3 estimates contrast_by_pre_wt(pupt_m, "beta")$emmeans
cat_pre_wt_log_scale = -2:
beta emmean SE df asymp.LCL asymp.UCL
-5 2.741 0.0941 Inf 2.557 2.93
0 2.790 0.0971 Inf 2.599 2.98
5 3.054 0.1088 Inf 2.841 3.27
cat_pre_wt_log_scale = 2:
beta emmean SE df asymp.LCL asymp.UCL
-5 1.061 0.0845 Inf 0.895 1.23
0 1.005 0.0805 Inf 0.847 1.16
5 0.866 0.0917 Inf 0.687 1.05
Results are averaged over the levels of: var_trt
Results are given on the log (not the response) scale.
Confidence level used: 0.95
$contrasts
cat_pre_wt_log_scale = -2:
contrast estimate SE df z.ratio p.value
(beta-5) - beta0 -0.0484 0.1018 Inf -0.476 0.8828
(beta-5) - beta5 -0.3132 0.1133 Inf -2.763 0.0158
beta0 - beta5 -0.2648 0.1136 Inf -2.331 0.0516
cat_pre_wt_log_scale = 2:
contrast estimate SE df z.ratio p.value
(beta-5) - beta0 0.0556 0.0875 Inf 0.636 0.8005
(beta-5) - beta5 0.1941 0.0915 Inf 2.122 0.0854
beta0 - beta5 0.1385 0.0928 Inf 1.492 0.2945
Results are averaged over the levels of: var_trt
Results are given on the log (not the response) scale.
P value adjustment: tukey method for comparing a family of 3 estimates 2.2.3 Misc stats
r2(pupt_m, tolerance = 10^-10)Random effect variances not available. Returned R2 does not account for random effects.# R2 for Mixed Models
Conditional R2: NA
Marginal R2: 0.8812.2.4 Check model
check_mod(pupt_m)
2.2.5 Plot
g3 <- marginal_effects(pupt_m, terms = c("cat_pre_wt_log_scale","beta")) %>%
ggplot(aes(x = unscalelog(d$cat_pre_wt_log)(cat_pre_wt_log_scale), y = yhat)) +
geom_ribbon(aes(ymax = upper, ymin = lower, fill = beta), alpha = 0.2) +
geom_line(aes(color = beta), linewidth = 2) +
geom_point(
data = pupt_m$frame,
aes(x = unscalelog(d$cat_pre_wt_log)(cat_pre_wt_log_scale), y = time_to_pupation, color = beta),
size = 3, alpha = 0.8
) +
theme_bw(base_size = 15) +
theme(legend.position = "top") +
scale_y_continuous(trans = "log2") +
scale_x_continuous(trans = "log10", labels = fancy_scientific) +
scale_color_brewer(type = "qual", aesthetics = c("fill","color")) +
labs(x = "Cat pre-weight (g)", y = "Time to pupation (days)",
color = expression(beta), fill = expression(beta))
g4 <- marginal_effects(pupt_m, terms = c("var_trt")) %>%
ggplot(aes(x = var_trt, y = yhat)) +
geom_point(
data = pupt_m$frame,
aes(x = var_trt, y = time_to_pupation, color = var_trt),
position = position_jitter(width = 0.2),
size = 3, alpha = 0.8
) +
geom_pointrange(
aes(ymax = upper, ymin = lower),
color = "black",
size = 1, linewidth = 2, shape = 1,
) +
theme_bw(base_size = 15) +
scale_y_continuous(trans = "log2") +
theme(legend.position = "top") +
scale_color_discrete(type = c("navy","skyblue"),
aesthetics = c("color"),
label = c("High", "Low")) +
labs(x = "Variation", y = "Time to pupation (days)",
color = "Variation")
ggarrange(g3, g4 + labs(y = ""))
2.3 Eclosure
2.3.1 Model summary
eclosure_m_full <- glmmTMB(
eclosed ~
(var_trt + beta) * cat_pre_wt_log_scale +
(1|session_id),
data = d %>%
filter(var_trt != "constant"),
family = binomial(link = "logit"),
); summary(eclosure_m_full) Family: binomial ( logit )
Formula:
eclosed ~ (var_trt + beta) * cat_pre_wt_log_scale + (1 | session_id)
Data: d %>% filter(var_trt != "constant")
AIC BIC logLik deviance df.resid
185.3 210.9 -83.6 167.3 119
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
session_id (Intercept) 0.1661 0.4076
Number of obs: 128, groups: session_id, 5
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 0.53530 0.41050 1.304 0.192
var_trtlow_var 0.01707 0.37507 0.046 0.964
beta0 -0.21544 0.46265 -0.466 0.641
beta5 -0.50749 0.45596 -1.113 0.266
cat_pre_wt_log_scale 0.10452 0.36089 0.290 0.772
var_trtlow_var:cat_pre_wt_log_scale 0.17526 0.36489 0.480 0.631
beta0:cat_pre_wt_log_scale 0.24778 0.43553 0.569 0.569
beta5:cat_pre_wt_log_scale 0.20159 0.44699 0.451 0.652car::Anova(eclosure_m_full, type = "III")Analysis of Deviance Table (Type III Wald chisquare tests)
Response: eclosed
Chisq Df Pr(>Chisq)
(Intercept) 1.7004 1 0.1922
var_trt 0.0021 1 0.9637
beta 1.2502 2 0.5352
cat_pre_wt_log_scale 0.0839 1 0.7721
var_trt:cat_pre_wt_log_scale 0.2307 1 0.6310
beta:cat_pre_wt_log_scale 0.3659 2 0.8328eclosure_m <- glmmTMB(
eclosed ~
(var_trt + beta) + cat_pre_wt_log_scale +
(1|session_id),
data = d %>%
filter(var_trt != "constant"),
family = binomial(link = "logit"),
); summary(eclosure_m) Family: binomial ( logit )
Formula:
eclosed ~ (var_trt + beta) + cat_pre_wt_log_scale + (1 | session_id)
Data: d %>% filter(var_trt != "constant")
AIC BIC logLik deviance df.resid
179.9 197.0 -84.0 167.9 122
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
session_id (Intercept) 0.1676 0.4094
Number of obs: 128, groups: session_id, 5
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 0.551286 0.413578 1.333 0.183
var_trtlow_var 0.005064 0.373402 0.014 0.989
beta0 -0.237692 0.460689 -0.516 0.606
beta5 -0.523521 0.456587 -1.147 0.252
cat_pre_wt_log_scale 0.330169 0.220049 1.500 0.134car::Anova(eclosure_m, type = "II")Analysis of Deviance Table (Type II Wald chisquare tests)
Response: eclosed
Chisq Df Pr(>Chisq)
var_trt 0.0002 1 0.9892
beta 1.3212 2 0.5165
cat_pre_wt_log_scale 2.2513 1 0.13352.3.2 Post-hoc contrasts
2.3.3 Misc stats
r2(eclosure_m, tolerance = 10^-10)# R2 for Mixed Models
Conditional R2: 0.092
Marginal R2: 0.0462.3.4 Check model
check_mod(eclosure_m)
2.3.5 Plot
2.4 Survival
2.4.1 Model summary
surv_m_full <- coxph(
Surv(surv_time, observed_dead) ~
(var_trt + beta) * cat_pre_wt_log_scale + I(cat_pre_wt_log_scale^2) +
strata(session_id),
data = d %>%
filter(var_trt != "constant") %>%
mutate(
beta = as.factor(as.character(beta))
),
); summary(surv_m_full)Call:
coxph(formula = Surv(surv_time, observed_dead) ~ (var_trt + beta) *
cat_pre_wt_log_scale + I(cat_pre_wt_log_scale^2) + strata(session_id),
data = d %>% filter(var_trt != "constant") %>% mutate(beta = as.factor(as.character(beta))))
n= 128, number of events= 105
coef exp(coef) se(coef) z Pr(>|z|)
var_trtlow_var -0.1115 0.8945 0.2124 -0.525 0.599727
beta0 -0.0896 0.9143 0.2579 -0.347 0.728306
beta5 -0.2521 0.7772 0.2605 -0.968 0.333144
cat_pre_wt_log_scale 0.9641 2.6224 0.2556 3.772 0.000162
I(cat_pre_wt_log_scale^2) 0.2450 1.2776 0.1447 1.693 0.090408
var_trtlow_var:cat_pre_wt_log_scale -0.1626 0.8499 0.1998 -0.814 0.415661
beta0:cat_pre_wt_log_scale -0.2230 0.8001 0.2309 -0.966 0.334053
beta5:cat_pre_wt_log_scale -0.4789 0.6195 0.2584 -1.854 0.063807
var_trtlow_var
beta0
beta5
cat_pre_wt_log_scale ***
I(cat_pre_wt_log_scale^2) .
var_trtlow_var:cat_pre_wt_log_scale
beta0:cat_pre_wt_log_scale
beta5:cat_pre_wt_log_scale .
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
exp(coef) exp(-coef) lower .95 upper .95
var_trtlow_var 0.8945 1.1179 0.5899 1.356
beta0 0.9143 1.0937 0.5515 1.516
beta5 0.7772 1.2867 0.4665 1.295
cat_pre_wt_log_scale 2.6224 0.3813 1.5892 4.328
I(cat_pre_wt_log_scale^2) 1.2776 0.7827 0.9621 1.696
var_trtlow_var:cat_pre_wt_log_scale 0.8499 1.1766 0.5746 1.257
beta0:cat_pre_wt_log_scale 0.8001 1.2499 0.5089 1.258
beta5:cat_pre_wt_log_scale 0.6195 1.6142 0.3734 1.028
Concordance= 0.637 (se = 0.047 )
Likelihood ratio test= 22 on 8 df, p=0.005
Wald test = 22.31 on 8 df, p=0.004
Score (logrank) test = 24.57 on 8 df, p=0.002drop1(surv_m_full, test = "Chisq")Single term deletions
Model:
Surv(surv_time, observed_dead) ~ (var_trt + beta) * cat_pre_wt_log_scale +
I(cat_pre_wt_log_scale^2) + strata(session_id)
Df AIC LRT Pr(>Chi)
<none> 459.50
I(cat_pre_wt_log_scale^2) 1 460.40 2.90 0.08841 .
strata(session_id) 0 787.87 328.37
var_trt:cat_pre_wt_log_scale 1 458.16 0.67 0.41395
beta:cat_pre_wt_log_scale 2 458.98 3.48 0.17532
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1surv_m <- coxph(
Surv(surv_time, observed_dead) ~
(var_trt + beta) + cat_pre_wt_log_scale + I(cat_pre_wt_log_scale^2) +
strata(session_id),
data = d %>%
filter(var_trt != "constant") %>%
mutate(
beta = as.factor(as.character(beta))
),
); summary(surv_m)Call:
coxph(formula = Surv(surv_time, observed_dead) ~ (var_trt + beta) +
cat_pre_wt_log_scale + I(cat_pre_wt_log_scale^2) + strata(session_id),
data = d %>% filter(var_trt != "constant") %>% mutate(beta = as.factor(as.character(beta))))
n= 128, number of events= 105
coef exp(coef) se(coef) z Pr(>|z|)
var_trtlow_var -0.07077 0.93167 0.20845 -0.340 0.734212
beta0 -0.07410 0.92858 0.25808 -0.287 0.774026
beta5 -0.17722 0.83759 0.25804 -0.687 0.492206
cat_pre_wt_log_scale 0.65947 1.93376 0.18164 3.631 0.000283 ***
I(cat_pre_wt_log_scale^2) 0.29584 1.34425 0.14124 2.095 0.036205 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
exp(coef) exp(-coef) lower .95 upper .95
var_trtlow_var 0.9317 1.0733 0.6192 1.402
beta0 0.9286 1.0769 0.5599 1.540
beta5 0.8376 1.1939 0.5051 1.389
cat_pre_wt_log_scale 1.9338 0.5171 1.3545 2.761
I(cat_pre_wt_log_scale^2) 1.3443 0.7439 1.0192 1.773
Concordance= 0.602 (se = 0.047 )
Likelihood ratio test= 18.08 on 5 df, p=0.003
Wald test = 18.29 on 5 df, p=0.003
Score (logrank) test = 19.42 on 5 df, p=0.002drop1(surv_m, test = "Chisq")Single term deletions
Model:
Surv(surv_time, observed_dead) ~ (var_trt + beta) + cat_pre_wt_log_scale +
I(cat_pre_wt_log_scale^2) + strata(session_id)
Df AIC LRT Pr(>Chi)
<none> 457.42
var_trt 1 455.53 0.12 0.7341150
beta 2 453.89 0.48 0.7883535
cat_pre_wt_log_scale 1 469.05 13.63 0.0002226 ***
I(cat_pre_wt_log_scale^2) 1 459.84 4.42 0.0355376 *
strata(session_id) 0 784.70 327.28
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 12.4.2 Post-hoc contrasts
2.4.3 Misc stats
r2(surv_m, tolerance = 10^-10) Nagelkerke's R2: 0.1352.4.4 Check model
cox.zph(surv_m) chisq df p
var_trt 1.60 1 0.21
beta 1.43 2 0.49
cat_pre_wt_log_scale 2.61 1 0.11
I(cat_pre_wt_log_scale^2) 1.05 1 0.31
GLOBAL 5.53 5 0.352.4.5 Plot
3 Structural Equation Model
3.1 Sub-models
3.1.1 RGR
subm_rgr <- glmmTMB(
RGR_scale ~
cat_pre_wt_log_scale +
cat_pre_wt_log_scale_sq +
sl_mean_obs_log_scale +
prop_explore_logit_scale *
var_toxic_12_scale +
mean_toxic_conc_scale +
area_herb_log_scale +
beta_numeric_scale +
(1|session_id),
data = d2
); summary(subm_rgr) Family: gaussian ( identity )
Formula: RGR_scale ~ cat_pre_wt_log_scale + cat_pre_wt_log_scale_sq +
sl_mean_obs_log_scale + prop_explore_logit_scale * var_toxic_12_scale +
mean_toxic_conc_scale + area_herb_log_scale + beta_numeric_scale +
(1 | session_id)
Data: d2
AIC BIC logLik deviance df.resid
119.5 149.0 -47.7 95.5 75
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
session_id (Intercept) 0.01245 0.1116
Residual 0.16744 0.4092
Number of obs: 87, groups: session_id, 5
Dispersion estimate for gaussian family (sigma^2): 0.167
Conditional model:
Estimate Std. Error z value
(Intercept) -0.02766 0.07257 -0.381
cat_pre_wt_log_scale -3.64523 0.41435 -8.797
cat_pre_wt_log_scale_sq -3.01946 0.45628 -6.618
sl_mean_obs_log_scale -0.13590 0.05626 -2.416
prop_explore_logit_scale 0.02370 0.04889 0.485
var_toxic_12_scale -0.01949 0.05115 -0.381
mean_toxic_conc_scale -0.25503 0.07147 -3.568
area_herb_log_scale 0.76522 0.09726 7.868
beta_numeric_scale -0.11265 0.04837 -2.329
prop_explore_logit_scale:var_toxic_12_scale 0.18687 0.04488 4.164
Pr(>|z|)
(Intercept) 0.703071
cat_pre_wt_log_scale < 2e-16 ***
cat_pre_wt_log_scale_sq 3.65e-11 ***
sl_mean_obs_log_scale 0.015710 *
prop_explore_logit_scale 0.627788
var_toxic_12_scale 0.703095
mean_toxic_conc_scale 0.000359 ***
area_herb_log_scale 3.61e-15 ***
beta_numeric_scale 0.019873 *
prop_explore_logit_scale:var_toxic_12_scale 3.13e-05 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 13.1.2 Variance in toxin (12 hour)
subm_var_toxic <- glmmTMB(
var_toxic_12_scale ~
sl_mean_obs_log_scale +
beta_numeric_scale +
var_high +
(1|session_id),
data = d2
); summary(subm_var_toxic) Family: gaussian ( identity )
Formula:
var_toxic_12_scale ~ sl_mean_obs_log_scale + beta_numeric_scale +
var_high + (1 | session_id)
Data: d2
AIC BIC logLik deviance df.resid
178.2 193.0 -83.1 166.2 81
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
session_id (Intercept) 6.229e-11 7.892e-06
Residual 3.956e-01 6.290e-01
Number of obs: 87, groups: session_id, 5
Dispersion estimate for gaussian family (sigma^2): 0.396
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -0.75075 0.09477 -7.922 2.34e-15 ***
sl_mean_obs_log_scale 0.21267 0.07751 2.744 0.00607 **
beta_numeric_scale -0.15689 0.06811 -2.303 0.02126 *
var_high 1.56021 0.13868 11.250 < 2e-16 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 13.1.3 Mean step length
subm_sl <- glmmTMB(
sl_mean_obs_log_scale ~
prop_explore_logit_scale +
cat_pre_wt_log_scale +
cat_pre_wt_log_scale_sq +
on_toxic_logit_scale +
(1|session_id),
data = d2,
); summary(subm_sl) Family: gaussian ( identity )
Formula:
sl_mean_obs_log_scale ~ prop_explore_logit_scale + cat_pre_wt_log_scale +
cat_pre_wt_log_scale_sq + on_toxic_logit_scale + (1 | session_id)
Data: d2
AIC BIC logLik deviance df.resid
223.4 240.7 -104.7 209.4 80
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
session_id (Intercept) 4.610e-10 2.147e-05
Residual 6.502e-01 8.064e-01
Number of obs: 87, groups: session_id, 5
Dispersion estimate for gaussian family (sigma^2): 0.65
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 0.06163 0.09646 0.639 0.52287
prop_explore_logit_scale -0.26175 0.08987 -2.912 0.00359 **
cat_pre_wt_log_scale 1.37948 0.70121 1.967 0.04915 *
cat_pre_wt_log_scale_sq 1.61095 0.75569 2.132 0.03303 *
on_toxic_logit_scale 0.21435 0.09489 2.259 0.02389 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 13.1.4 Mean toxin ingested
subm_toxin_ingested <- glmmTMB(
mean_toxic_conc_scale ~
beta_numeric_scale +
on_toxic_logit_scale +
var_high +
(1|session_id),
data = d2
); summary(subm_toxin_ingested) Family: gaussian ( identity )
Formula:
mean_toxic_conc_scale ~ beta_numeric_scale + on_toxic_logit_scale +
var_high + (1 | session_id)
Data: d2
AIC BIC logLik deviance df.resid
145.7 160.5 -66.8 133.7 81
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
session_id (Intercept) 6.461e-10 2.542e-05
Residual 2.722e-01 5.217e-01
Number of obs: 87, groups: session_id, 5
Dispersion estimate for gaussian family (sigma^2): 0.272
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 0.18349 0.07861 2.334 0.0196 *
beta_numeric_scale -0.11683 0.05808 -2.012 0.0443 *
on_toxic_logit_scale 0.28426 0.06393 4.447 8.72e-06 ***
var_high -0.80646 0.11613 -6.945 3.79e-12 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 13.1.5 Mean time on toxin
subm_on_toxic <- glmmTMB(
on_toxic_logit_scale ~
var_high + beta_numeric_scale * cat_pre_wt_log_scale +
(1|session_id),
family = gaussian(),
data = d2
); summary(subm_on_toxic) Family: gaussian ( identity )
Formula:
on_toxic_logit_scale ~ var_high + beta_numeric_scale * cat_pre_wt_log_scale +
(1 | session_id)
Data: d2
AIC BIC logLik deviance df.resid
230.0 247.3 -108.0 216.0 80
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
session_id (Intercept) 0.01667 0.1291
Residual 0.68757 0.8292
Number of obs: 87, groups: session_id, 5
Dispersion estimate for gaussian family (sigma^2): 0.688
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 0.22509 0.14140 1.592 0.1114
var_high -0.41809 0.18032 -2.319 0.0204 *
beta_numeric_scale -0.12409 0.09766 -1.271 0.2039
cat_pre_wt_log_scale -0.19306 0.13787 -1.400 0.1614
beta_numeric_scale:cat_pre_wt_log_scale -0.23850 0.11732 -2.033 0.0421 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 13.1.6 Area consumption
subm_herb <- glmmTMB(
area_herb_log_scale ~
var_high * cat_pre_wt_log_scale +
(1|session_id),
data = d2,
); summary(subm_herb) Family: gaussian ( identity )
Formula: area_herb_log_scale ~ var_high * cat_pre_wt_log_scale + (1 |
session_id)
Data: d2
AIC BIC logLik deviance df.resid
117.9 132.7 -53.0 105.9 81
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
session_id (Intercept) 3.501e-11 5.917e-06
Residual 1.978e-01 4.448e-01
Number of obs: 87, groups: session_id, 5
Dispersion estimate for gaussian family (sigma^2): 0.198
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) 0.04339 0.06963 0.623 0.53313
var_high 0.27072 0.10127 2.673 0.00751 **
cat_pre_wt_log_scale 0.55603 0.08378 6.637 3.2e-11 ***
var_high:cat_pre_wt_log_scale -0.26289 0.11558 -2.275 0.02293 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 13.1.7 Exploration
subm_exp <- glmmTMB(
prop_explore_logit_scale ~
var_high + beta_numeric_scale * cat_pre_wt_log_scale +
(1|session_id),
family = gaussian(),
data = d2
); summary(subm_exp) Family: gaussian ( identity )
Formula:
prop_explore_logit_scale ~ var_high + beta_numeric_scale * cat_pre_wt_log_scale +
(1 | session_id)
Data: d2
AIC BIC logLik deviance df.resid
250.1 267.3 -118.0 236.1 80
Random effects:
Conditional model:
Groups Name Variance Std.Dev.
session_id (Intercept) 3.781e-10 1.944e-05
Residual 8.830e-01 9.397e-01
Number of obs: 87, groups: session_id, 5
Dispersion estimate for gaussian family (sigma^2): 0.883
Conditional model:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -0.17231 0.14397 -1.197 0.2314
var_high 0.15661 0.20355 0.769 0.4416
beta_numeric_scale -0.14059 0.11040 -1.273 0.2029
cat_pre_wt_log_scale -0.04358 0.12362 -0.353 0.7244
beta_numeric_scale:cat_pre_wt_log_scale 0.25977 0.13242 1.962 0.0498 *
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 13.2 d-separation test
sem_fit <- psem(
subm_rgr,
subm_herb,
subm_var_toxic,
subm_on_toxic,
subm_toxin_ingested,
subm_exp,
subm_sl,
var_toxic_12_scale %~~% on_toxic_logit_scale, # correlated errors
data = d2
)
sem_summary <- suppressWarnings(suppressMessages(
summary_psem(sem_fit,
no_standardize_x = c(),
.progressBar = FALSE,
direction = "var_toxic_12_scale -> area_herb_log_scale")
))Random effect variances not available. Returned R2 does not account for random effects.
Random effect variances not available. Returned R2 does not account for random effects.
Random effect variances not available. Returned R2 does not account for random effects.
Random effect variances not available. Returned R2 does not account for random effects.
Random effect variances not available. Returned R2 does not account for random effects.
Random effect variances not available. Returned R2 does not account for random effects.
Random effect variances not available. Returned R2 does not account for random effects.
Random effect variances not available. Returned R2 does not account for random effects.
Random effect variances not available. Returned R2 does not account for random effects.
Random effect variances not available. Returned R2 does not account for random effects.sem_summary$Cstat Fisher.C df P.Value
1 42.469 42 0.451knitr::kable(sem_summary$dTable)| Independ.Claim | Test.Type | DF | Crit.Value | P.Value | |
|---|---|---|---|---|---|
| mean_toxic_conc_scale ~ cat_pre_wt_log_scale + … | coef | 87 | -1.5712 | 0.1161 | |
| var_toxic_12_scale ~ cat_pre_wt_log_scale + … | coef | 87 | -0.1892 | 0.8499 | |
| prop_explore_logit_scale ~ cat_pre_wt_log_scale_sq + … | coef | 87 | 0.8624 | 0.3884 | |
| mean_toxic_conc_scale ~ cat_pre_wt_log_scale_sq + … | coef | 87 | 1.8002 | 0.0718 | |
| area_herb_log_scale ~ cat_pre_wt_log_scale_sq + … | coef | 87 | -0.7057 | 0.4804 | |
| on_toxic_logit_scale ~ cat_pre_wt_log_scale_sq + … | coef | 87 | 0.7944 | 0.4269 | |
| var_toxic_12_scale ~ cat_pre_wt_log_scale_sq + … | coef | 87 | 0.3081 | 0.7580 | |
| sl_mean_obs_log_scale ~ beta_numeric_scale + … | coef | 87 | 1.1020 | 0.2704 | |
| area_herb_log_scale ~ beta_numeric_scale + … | coef | 87 | 0.0200 | 0.9840 | |
| RGR_scale ~ var_high + … | coef | 87 | 0.6991 | 0.4845 | |
| sl_mean_obs_log_scale ~ var_high + … | coef | 87 | -1.5883 | 0.1122 | |
| on_toxic_logit_scale ~ RGR_scale + … | coef | 87 | -1.3664 | 0.1718 | |
| mean_toxic_conc_scale ~ sl_mean_obs_log_scale + … | coef | 87 | 0.7459 | 0.4557 | |
| area_herb_log_scale ~ sl_mean_obs_log_scale + … | coef | 87 | 0.8564 | 0.3918 | |
| mean_toxic_conc_scale ~ prop_explore_logit_scale + … | coef | 87 | 1.0591 | 0.2895 | |
| area_herb_log_scale ~ prop_explore_logit_scale + … | coef | 87 | 0.3920 | 0.6951 | |
| on_toxic_logit_scale ~ prop_explore_logit_scale + … | coef | 87 | 0.2361 | 0.8134 | |
| var_toxic_12_scale ~ prop_explore_logit_scale + … | coef | 87 | 1.1110 | 0.2666 | |
| area_herb_log_scale ~ mean_toxic_conc_scale + … | coef | 87 | -0.2457 | 0.8059 | |
| var_toxic_12_scale ~ mean_toxic_conc_scale + … | coef | 87 | -1.3086 | 0.1907 | |
| on_toxic_logit_scale ~ area_herb_log_scale + … | coef | 87 | -0.7213 | 0.4707 |
3.3 Model summary
3.3.1 Model stats
sem_summary$ChiSqNULLsem_summary$AIC AIC AICc K n
1 1264.813 1276.432 51 87sem_summary$R2 Response family link method Marginal Conditional
1 RGR_scale gaussian identity none 0.75 0.77
2 area_herb_log_scale gaussian identity none 0.42 NA
3 var_toxic_12_scale gaussian identity none 0.60 NA
4 on_toxic_logit_scale gaussian identity none 0.19 0.21
5 mean_toxic_conc_scale gaussian identity none 0.55 NA
6 prop_explore_logit_scale gaussian identity none 0.06 NA
7 sl_mean_obs_log_scale gaussian identity none 0.19 NA3.3.2 Standardized coefficients
knitr::kable(
cbind(sem_summary$coefficients,"arrow_size"=(abs(as.numeric(sem_summary$coefficients$Std.Estimate))*20))
)| Response | Predictor | Estimate | Std.Error | DF | Crit.Value | P.Value | Std.Estimate | star | arrow_size |
|---|---|---|---|---|---|---|---|---|---|
| RGR_scale | cat_pre_wt_log_scale | -3.6452 | 0.4144 | 87 | -8.7974 | 0.0000 | -3.5981 | *** | 71.962 |
| RGR_scale | cat_pre_wt_log_scale_sq | -3.0195 | 0.4563 | 87 | -6.6176 | 0.0000 | -2.7708 | *** | 55.416 |
| RGR_scale | sl_mean_obs_log_scale | -0.1359 | 0.0563 | 87 | -2.4156 | 0.0157 | -0.1448 | * | 2.896 |
| RGR_scale | prop_explore_logit_scale | 0.0237 | 0.0489 | 87 | 0.4848 | 0.6278 | 0.0274 | 0.548 | |
| RGR_scale | var_toxic_12_scale | -0.0195 | 0.0511 | 87 | -0.3811 | 0.7031 | -0.0230 | 0.460 | |
| RGR_scale | mean_toxic_conc_scale | -0.2550 | 0.0715 | 87 | -3.5683 | 0.0004 | -0.2353 | *** | 4.706 |
| RGR_scale | area_herb_log_scale | 0.7652 | 0.0973 | 87 | 7.8679 | 0.0000 | 0.5338 | *** | 10.676 |
| RGR_scale | beta_numeric_scale | -0.1126 | 0.0484 | 87 | -2.3287 | 0.0199 | -0.1337 | * | 2.674 |
| RGR_scale | prop_explore_logit_scale:var_toxic_12_scale | 0.1869 | 0.0449 | 87 | 4.1640 | 0.0000 | 0.2280 | *** | 4.560 |
| area_herb_log_scale | var_high | 0.2707 | 0.1013 | 87 | 2.6733 | 0.0075 | 0.2315 | ** | 4.630 |
| area_herb_log_scale | cat_pre_wt_log_scale | 0.5560 | 0.0838 | 87 | 6.6371 | 0.0000 | 0.7868 | *** | 15.736 |
| area_herb_log_scale | var_high:cat_pre_wt_log_scale | -0.2629 | 0.1156 | 87 | -2.2745 | 0.0229 | -0.2795 | * | 5.590 |
| var_toxic_12_scale | sl_mean_obs_log_scale | 0.2127 | 0.0775 | 87 | 2.7437 | 0.0061 | 0.1918 | ** | 3.836 |
| var_toxic_12_scale | beta_numeric_scale | -0.1569 | 0.0681 | 87 | -2.3033 | 0.0213 | -0.1575 | * | 3.150 |
| var_toxic_12_scale | var_high | 1.5602 | 0.1387 | 87 | 11.2503 | 0.0000 | 0.7875 | *** | 15.750 |
| on_toxic_logit_scale | var_high | -0.4181 | 0.1803 | 87 | -2.3186 | 0.0204 | -0.2235 | * | 4.470 |
| on_toxic_logit_scale | beta_numeric_scale | -0.1241 | 0.0977 | 87 | -1.2706 | 0.2039 | -0.1320 | 2.640 | |
| on_toxic_logit_scale | cat_pre_wt_log_scale | -0.1931 | 0.1379 | 87 | -1.4003 | 0.1614 | -0.1708 | 3.416 | |
| on_toxic_logit_scale | beta_numeric_scale:cat_pre_wt_log_scale | -0.2385 | 0.1173 | 87 | -2.0329 | 0.0421 | -0.2116 | * | 4.232 |
| mean_toxic_conc_scale | beta_numeric_scale | -0.1168 | 0.0581 | 87 | -2.0116 | 0.0443 | -0.1503 | * | 3.006 |
| mean_toxic_conc_scale | on_toxic_logit_scale | 0.2843 | 0.0639 | 87 | 4.4467 | 0.0000 | 0.3438 | *** | 6.876 |
| mean_toxic_conc_scale | var_high | -0.8065 | 0.1161 | 87 | -6.9447 | 0.0000 | -0.5214 | *** | 10.428 |
| prop_explore_logit_scale | var_high | 0.1566 | 0.2035 | 87 | 0.7694 | 0.4416 | 0.0808 | 1.616 | |
| prop_explore_logit_scale | beta_numeric_scale | -0.1406 | 0.1104 | 87 | -1.2734 | 0.2029 | -0.1443 | 2.886 | |
| prop_explore_logit_scale | cat_pre_wt_log_scale | -0.0436 | 0.1236 | 87 | -0.3526 | 0.7244 | -0.0372 | 0.744 | |
| prop_explore_logit_scale | beta_numeric_scale:cat_pre_wt_log_scale | 0.2598 | 0.1324 | 87 | 1.9617 | 0.0498 | 0.2224 | * | 4.448 |
| sl_mean_obs_log_scale | prop_explore_logit_scale | -0.2618 | 0.0899 | 87 | -2.9124 | 0.0036 | -0.2840 | ** | 5.680 |
| sl_mean_obs_log_scale | cat_pre_wt_log_scale | 1.3795 | 0.7012 | 87 | 1.9673 | 0.0491 | 1.2778 | * | 25.556 |
| sl_mean_obs_log_scale | cat_pre_wt_log_scale_sq | 1.6109 | 0.7557 | 87 | 2.1318 | 0.0330 | 1.3873 | * | 27.746 |
| sl_mean_obs_log_scale | on_toxic_logit_scale | 0.2144 | 0.0949 | 87 | 2.2589 | 0.0239 | 0.2244 | * | 4.488 |
| ~~var_toxic_12_scale | ~~on_toxic_logit_scale | 0.5323 | - | 87 | 5.7626 | 0.0000 | 0.5323 | *** | 10.646 |
3.3.3 Post hoc contrast
contrast_by_pre_wt(subm_herb, "var_high", "trt")$emmeans
cat_pre_wt_log_scale = -2:
var_high emmean SE df lower.CL upper.CL
0 -1.069 0.200 81 -1.467 -0.670
1 -0.272 0.198 81 -0.666 0.122
cat_pre_wt_log_scale = 2:
var_high emmean SE df lower.CL upper.CL
0 1.155 0.161 81 0.836 1.475
1 0.900 0.150 81 0.603 1.198
Confidence level used: 0.95
$contrasts
cat_pre_wt_log_scale = -2:
contrast estimate SE df t.ratio p.value
var_high1 - var_high0 0.797 0.281 81 2.830 0.0059
cat_pre_wt_log_scale = 2:
contrast estimate SE df t.ratio p.value
var_high1 - var_high0 -0.255 0.219 81 -1.162 0.2485sd(subm_herb$frame$var_high)/sd(subm_herb$frame$area_herb_log_scale)[1] 0.8549509emmeans::emtrends(subm_exp,
~ beta_numeric_scale | cat_pre_wt_log_scale,
var = "beta_numeric_scale",
at = list(cat_pre_wt_log_scale = c(-2,2)))cat_pre_wt_log_scale = -2:
beta_numeric_scale beta_numeric_scale.trend SE df lower.CL upper.CL
7.88e-17 -0.660 0.323 80 -1.30 -0.0173
cat_pre_wt_log_scale = 2:
beta_numeric_scale beta_numeric_scale.trend SE df lower.CL upper.CL
7.88e-17 0.379 0.246 80 -0.11 0.8677
Results are averaged over the levels of: var_high
Confidence level used: 0.95 sd(subm_exp$frame$beta_numeric_scale) * sd(subm_exp$frame$prop_explore_logit_scale)[1] 0.974435emmeans::emtrends(subm_on_toxic,
~ beta_numeric_scale | cat_pre_wt_log_scale,
var = "beta_numeric_scale",
at = list(cat_pre_wt_log_scale = c(-2,2)))cat_pre_wt_log_scale = -2:
beta_numeric_scale beta_numeric_scale.trend SE df lower.CL upper.CL
7.88e-17 0.353 0.286 80 -0.217 0.923
cat_pre_wt_log_scale = 2:
beta_numeric_scale beta_numeric_scale.trend SE df lower.CL upper.CL
7.88e-17 -0.601 0.217 80 -1.034 -0.169
Results are averaged over the levels of: var_high
Confidence level used: 0.95 sd(subm_on_toxic$frame$beta_numeric_scale) * sd(subm_on_toxic$frame$on_toxic_logit_scale)[1] 0.94033613.4 SEM plots
knitr::include_graphics(paste0(getwd(), "/graphs/manuscript1_figures/SEM_simplified.png"), rel_path = FALSE)
4 Session Info
sessionInfo()R version 4.3.1 (2023-06-16 ucrt)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19045)
Matrix products: default
locale:
[1] LC_COLLATE=English_United States.utf8
[2] LC_CTYPE=English_United States.utf8
[3] LC_MONETARY=English_United States.utf8
[4] LC_NUMERIC=C
[5] LC_TIME=English_United States.utf8
time zone: America/New_York
tzcode source: internal
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] spat1f_0.0.1 herbivar_0.2.1 forcats_0.5.1 stringr_1.5.0
[5] dplyr_1.1.2 purrr_1.0.1 readr_2.1.2 tidyr_1.3.0
[9] tibble_3.2.1 tidyverse_1.3.1 extraDistr_1.9.1 imager_0.42.13
[13] magrittr_2.0.3 doSNOW_1.0.20 snow_0.4-4 iterators_1.0.14
[17] foreach_1.5.2 glmmTMB_1.1.7 tidybayes_3.0.2 performance_0.10.4
[21] sjPlot_2.8.15 ggpubr_0.6.0 ggplot2_3.4.2 piecewiseSEM_2.3.0
[25] survival_3.3-1
loaded via a namespace (and not attached):
[1] svUnit_1.0.6 splines_4.3.1 later_1.3.0
[4] bitops_1.0-7 cellranger_1.1.0 polyclip_1.10-0
[7] reprex_2.0.1 lifecycle_1.0.3 sf_1.0-8
[10] Rdpack_2.3.1 rstatix_0.7.2 doParallel_1.0.17
[13] rprojroot_2.0.3 vroom_1.5.7 processx_3.6.1
[16] lattice_0.21-8 MASS_7.3-60 insight_0.19.3
[19] readbitmap_0.1.5 ggdist_3.1.1 backports_1.4.1
[22] sass_0.4.1 rmarkdown_2.24 jquerylib_0.1.4
[25] yaml_2.3.5 remotes_2.4.2 httpuv_1.6.5
[28] gap_1.2.3-6 sp_1.5-0 sessioninfo_1.2.2
[31] pkgbuild_1.3.1 cowplot_1.1.1 DBI_1.1.3
[34] minqa_1.2.4 RColorBrewer_1.1-3 lubridate_1.8.0
[37] DHARMa_0.4.6 multcomp_1.4-25 abind_1.4-5
[40] pkgload_1.3.2 rvest_1.0.2 imagerExtra_2.0.0
[43] TH.data_1.1-1 tensorA_0.36.2 sandwich_3.0-2
[46] spatstat.utils_3.0-1 units_0.8-0 fitdistrplus_1.1-11
[49] codetools_0.2-19 MuMIn_1.47.5 xml2_1.3.3
[52] tidyselect_1.2.0 ggeffects_1.2.3 farver_2.1.0
[55] lme4_1.1-29 matrixStats_1.2.0 stats4_4.3.1
[58] jsonlite_1.8.8 e1071_1.7-11 ellipsis_0.3.2
[61] bmp_0.3 emmeans_1.7.5 tools_4.3.1
[64] Rcpp_1.0.8.3 glue_1.6.2 mgcv_1.8-40
[67] xfun_0.39 distributional_0.3.0 usethis_2.1.6
[70] withr_2.5.0 numDeriv_2016.8-1.1 fastmap_1.1.0
[73] boot_1.3-28.1 fansi_1.0.3 callr_3.7.0
[76] digest_0.6.29 R6_2.5.1 mime_0.12
[79] estimability_1.3 colorspace_2.0-3 spatstat.data_3.0-0
[82] jpeg_0.1-9 see_0.7.1 DiagrammeR_1.0.9
[85] utf8_1.2.2 generics_0.1.2 class_7.3-22
[88] prettyunits_1.1.1 httr_1.4.3 htmlwidgets_1.5.4
[91] pkgconfig_2.0.3 gtable_0.3.0 circular_0.5-0
[94] moveHMM_1.9 htmltools_0.5.2 carData_3.0-5
[97] profvis_0.3.7 fftwtools_0.9-11 TMB_1.9.4
[100] scales_1.2.1 png_0.1-7 posterior_1.2.2
[103] CircStats_0.2-6 knitr_1.43 rstudioapi_0.13
[106] geosphere_1.5-18 tzdb_0.3.0 coda_0.19-4
[109] visNetwork_2.1.0 checkmate_2.1.0 nlme_3.1-162
[112] nloptr_2.0.3 proxy_0.4-27 cachem_1.0.6
[115] zoo_1.8-12 sjlabelled_1.2.0 KernSmooth_2.23-20
[118] parallel_4.3.1 miniUI_0.1.1.1 desc_1.4.1
[121] pillar_1.9.0 grid_4.3.1 vctrs_0.6.3
[124] urlchecker_1.0.1 promises_1.2.0.1 car_3.1-2
[127] dbplyr_2.2.0 arrayhelpers_1.1-0 xtable_1.8-4
[130] evaluate_0.15 amt_0.2.1.0 mvtnorm_1.1-3
[133] cli_3.6.1 compiler_4.3.1 rlang_1.1.1
[136] crayon_1.5.1 ggsignif_0.6.3 labeling_0.4.2
[139] ggmap_3.0.2 modelr_0.1.8 classInt_0.4-7
[142] RcppArmadillo_0.11.2.0.0 ps_1.7.1 plyr_1.8.7
[145] sjmisc_2.8.9 fs_1.5.2 gap.datasets_0.0.5
[148] stringi_1.7.12 deldir_1.0-6 assertthat_0.2.1
[151] munsell_0.5.0 tiff_0.1-11 spatstat.geom_3.0-3
[154] devtools_2.4.5 bayestestR_0.13.1 Matrix_1.5-4.1
[157] sjstats_0.18.2 hms_1.1.1 bit64_4.0.5
[160] qgam_1.3.4 RgoogleMaps_1.4.5.3 shiny_1.7.1
[163] highr_0.9 haven_2.5.0 rbibutils_2.2.8
[166] igraph_1.3.2 broom_1.0.5 memoise_2.0.1
[169] bslib_0.3.1 bit_4.0.4 readxl_1.4.0